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Abstract 

Electrical transport in semiconductor superlattices is studied within a fully 
self-consistent quantum transport model based on nonequilibrium Green func- 
tions, including phonon and impurity scattering. We compute both the drift 
velocity-field relation and the momentum distribution function covering the 
whole field range from linear response to negative differential conductivity. 
The quantum results are compared with the respective results obtained from 
a Monte Carlo solution of the Boltzmann equation. Our analysis thus sets the 
limits of validity for the semiclassical theory in a nonlinear transport situation 
in the presence of inelastic scattering. 
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Quantum mechanical description of electrical transport in strong electric fields is a noto- 
riously difficult subject. As the distribution function of the electrons deviates strongly from 
equilibrium, the standard approaches of linear response do not apply. In some circumstances 
heated distribution functions may be useful, but in principle such an assumption should be 
justified by an underlying theory. Usually this problem of hot electrons is treated within 
the semiclassical Boltzmann transport equation (BTE) which can be solved to a desired 
degree of numerical accuracy by Monte-Carlo simulations (MC) pj. This way both the 
distribution function and the current density can be obtained for arbitrary field strengths. 
Nevertheless one has to be aware of the severe assumptions implied by the use of the BTE. 
The electrons are viewed as classical particles with a dispersion relation given by the band 
structure and scattering is treated using Fermi's golden rule where energy conservation is 
strict. To recover the BTE from a quantum transport theory, such as nonequilibrium Green 
functions (NGF) [|], several assumptions are required (see e.g. H||): (i) spectral functions 
with finite width are replaced by ^-functions for free particles, (ii) the nonequilibrium Green 
function is assumed to be expressible in terms of the momentum distribution function, (iii) 
retardation effects are neglected to reproduce the Markovian collision integral ||. Attempts 
to relax some of these assumptions have been made in earlier studies [||,[7| 5 and quantum 
corrections to distribution functions have been reported. Nevertheless we are not aware of 
any direct comparison of the results from BTE with a full quantum transport theory far 
from equilibrium, which is the task of this paper. 

Semiconductor superlattices (SL) || provide a unique opportunity to study effects related 
to quantum transport because the width of the miniband can be tailored by the choice of 
the barrier and well widths as well as the material composition. For sufficiently large electric 
fields negative differential conductivity (NDC) appears which can be understood within the 
semiclassical theory by electrons traversing the whole miniband || thus performing Bloch- 
oscillations both in momentum and real space. This gives a distribution function which 
is far from any kind of thermal equilibrium. Many papers have analyzed the solution of 
the BTE in this miniband transport regime ||10|| . Alternatively, in the NDC region the 



electrical transport can be formulated in terms of Wannier-Stark hopping [11 . Finally, for 



small miniband width A and strong scattering, the transport can be viewed as a sequential 



tunneling process between adjacent wells [12 



While these simplified approaches have proven to be useful to analyze different experi- 
mental situations, it is clear that a complete description of transport in SLs requires a more 
sophisticated quantum mechanical treatment such as density matrix theory [0,0] or NGF 



fl5| . Recently it was shown that a calculation based on NGF reproduces the results of 
the simplified approaches and determines their respective ranges of validity [IE|. All these 
approaches employed a specific (heated) thermal distribution to model the in-scattering 
processes. While this can be a reasonable approximation for the evaluation of averaged 
quantities such as the current density, such an assumption is certainly not justified if the 
nonequilibrium electron distribution function itself is of interest. 

In this work we calculate the drift velocity-field relation as well as the electron dis- 
tribution function both in a quantum transport model based on NGF and within the BTE 
under stationary conditions. Identical system parameters and scattering matrix elements for 
impurity and phonon scattering are used. Our NGF calculation provides us with a full self- 
consistent solution of the transport problem within the self-consistent Born-approximation 



for the scattering. (Thus, effects of higher order scattering such as weak localization are 

neglected.) This approach is similar to recent calculations for the resonant tunneling diode 

T7| . The BTE is solved by a MC-simulation. Its results can be easily interpreted in terms of 



standard concepts such as electron heating and Bloch-oscillating electrons. We find excellent 
agreement between the drift velocities calculated semiclassically and quantum mechanically 



in the range of validity for the miniband transport discussed in Refs. ||16| , p^8| . On the other 
hand, the failure of the BTE becomes obvious if the potential drop per period eFd or the 
scattering width h/r scautt become larger than the miniband width A. 

We consider a semiconductor superlattice with period d and restrict ourselves to nearest 
neighbor coupling. The band structure reads E(q, k) = 27\ cos(qd) + E^ with T\ = —A/4 
and Ek = h 2 k 2 /2m. Here q denotes the Bloch vector in the growth direction, which is 
restricted to —n/d < q < rc/d, and k is the two-dimensional Bloch vector perpendicular to 
the growth direction, m = 0.067m e is the effective mass of the conduction band of GaAs. We 
restrict ourselves to the lowest miniband neglecting all phenomena related to intersubband 
processes (e.g., the second peak in the drift velocity- field relation at higher fields ||[12|)- 

Three types of scattering processes are included: Impurity scattering at 5-potentials with 
density Nj (per unit area) and constant matrix element Vj mp leading to a scattering rate 
l/'Hmp = NdTrV^pPo/h, where p is the two-dimensional density of states. (The impact of 
further elastic scattering processes like interface roughness can be taken into account by 
applying a reduced value of r imp .) Optical phonons of energy ftw = 36 meV with a constant 
matrix element M Q . We choose M such that the rate for spontaneous emission of a phonon 
(if allowed) is given by l/r = nM 2 p /h = 8 ps _1 . These values are realistic for GaAs. In 
order to achieve energy relaxation for energies lower than hw we mimic acoustic phonons 
by a second phonon with constant energy frw^. The ratio uj /u a should be irrational to avoid 
spurious resonances and we choose fojj & = ftw (\^E — 1)/10 ~ 4.4498 meV. The constant 



matrix element M a was chosen to yield a scattering rate l/r a = 200 ns [fL9 l- These 
matrix elements are used in the scattering term of the BTE assuming a thermal occupation 
-Wo/a — [exp(/iu; 0/ / a /£;BT) — 1] _1 for the phonon modes with lattice temperature T. We solve 
the BTE by a MC procedure and obtain the semiclassical distribution function /sc (<?, k) 
as well as the average drift velocity thrift as a function of the electric field F applied to the 
SL. 



Following Ref. [jnj we use a basis of Wannier functions ^(z — nd) (localized in well n) 
for the NGF calculation. The task is to evaluate the retarded and lesser Green function 
G^ n (t,t',k) = -*6(t-t')({a m (t,k),c4(t',k)}> and G<Jt,t',k) = i<4(f,k)a m (t,k)), re- 
spectively. Here a^(t,k) and a n (t,k) are the creation and annihilation operators for the 
state ^(z — nd)e l ( kr ^ /A and {A, B} denotes the anticommutator. Using the Dyson and 



Keldysh equation [Eqs. (13,15) of Ref. fig] ] the Green functions can be calculated for given 
self-energies S^ et and £^ (which are diagonal in the well index as short range scattering po- 
tentials are assumed). In Ref. |]T6| only elastic scattering was considered and an equilibrium 
approximation for E^ was made to ensure energy relaxation. Here, instead, we calculate 
both S^ et and E^ self-consistently. For impurity scattering we use 



k' 
and for phonon scattering (see, e.g., Ch. 4.3 of Ref. 



^n.impl^) _ a Z^ Mmp^n,!! \y J k ) (1) 



\M I 2 
K&) = HrT, {NoG< n {£ - ^ , k') + (N + l)G< n (£ + hv Q , k')} (2) 



k' 



\M r r 
K e l(£) = ^p £ {(N + 1)G%J£ - t*> , k') + iVoG^ B (£ + &*, k') 



£' - hu + i0+ £> + hw Q + iO- 



(3) 



In our numerical calculation we ignore the real part of the last term containing the energy 
integral in E^*, which renormalizes the energy slightly. The contribution due to 'acoustic' 
phonons u> a is treated analogously. Note that the constant matrix elements lead to k- 
independent self-energies which implies a significant reduction in the computational effort. 
In our calculation we start with a guess for the self-energies, calculate the Green-functions 
and obtain a new set of self-energies from Eqs. (|l]-0). This procedure is continued iteratively 
until the new self-energies deviate by less than 0.5% from the preceding ones. We checked 
the convergence by comparing the currents starting from different initial guesses and found 
that even a stricter bound had to be used for low electric fields. Finally, we calculate the 
quantum momentum distribution function 

f QM (q, k) = JL J d S £ e- ih « d G< (£, k) (4) 

from the nondiagonal elements of G < and the current density 

'^p/Ow'- 1 '^ (5) 



which is equivalent to Eq. (10) of Ref. [J16 



Results: We first consider a wide-band superlattice with A = 20.3 meV and period 
d — 5.1 nm. We assume an average carrier density of 10 16 cm -3 and T = 77 K. The impurity 
scattering rate is taken to be l/ri mp = 3ps _1 . Fig. [1] shows the calculated drift- velocity 
versus electric field. We find that the characteristics calculated from NGF and BTE are in 
excellent agreement for eFd ;$ A/2 as expected from the analysis of Ref. [IB| because h/r ~ 



2.3meV <^ A/2. The shape of the Vdntt(F) relation significantly deviates from the simple 
Esaki-Tsu shape v drift oc F/(F 2 +F^) with F c = h/edr. A linear part is only observed for very 
low electric fields. Here the distribution function (see Fig. 0(a)) can be viewed as a distorted 
thermal equilibrium function. Thus the standard theory of linear response || makes sense 
yielding a linear part of the characteristics (weak localization effects may affect the good 
agreement between quantum transport and BTE at low temperatures). In Fig. |3](a) we 
have shown the respective distribution versus k = |k|, where f(k) = d/(2ir) Jj 7r , d dqf(q,\s). 
We find fsc{k) ~ Cexp(— E^/hsT) from the BTE, where C is a normalization constant. 
In the NGF calculation this behavior is only seen in the range Ef. < 20 meV. For higher 
values of E^ the quantum mechanical result is larger than the semiclassical one as energy 
broadening leads to the power-law /qm(&) ~ CksTT /(2tiED in the momentum distribution 
function, where T is the total scattering width. 

If the electric field increases, electron heating becomes important. For moderate fields 
the distribution function resembles a distorted equilibrium f eq (q, k) oc exp[— E(q, k)/A; B T e ] 



with an increased electron temperature T e p» 140 K for eFd = 0.3 meV, see Fig. 0(b). 
This suppresses the mobility and causes a sublinear increase of the current. Close to the 
maximum at eFd = 2 meV the distribution function strongly deviates from any kind of 
equilibrium in g-space (see Fig. |](c)), but the fc-dependence can be still viewed as a heated 
distribution (see Fig. 0(b)) with T e « 190 K for Ef. < fko . The results look similar in 
the NDC region for eFd = 10 meV (not shown here). In all cases (with eFd < A/2) the 
distribution functions from BTE agree well with the result from a full quantum mechanical 
description. Thus, the BTE can be viewed adequate in this parameter range. 

The situation chances dramatically for eFd > A, see Figs. |2|(d),|3](c). As the electrons 
can perform several Bloch-oscillations in the semiclassical picture, the distribution function 
is almost flat within the Brillouin zone of the miniband. The latter holds for the NGF 
result as well. However, the absolute values of the distribution functions differ significantly. 
The reason is the modification in scattering processes due to the presence of the electric 
field, leading to significant deviations in the distribution function. Therefore, also the drift 
velocities deviate significantly in this field region, see Fig. [I]. Finally, note the phonon- 
resonance Jl^^Oj at eFd = huj Q in the NGF result for the velocity-field relation. This feature 
cannot be recovered from the BTE where the field does not appear as an energy scale in the 
scattering term. The strong change of /qm(^) close to the phonon resonance in the NGF 
calculation is shown in Fig. 0(d). In the high-field regions one typically encounters patterns 
on the energy scales fuv and eFd (due to the formation of the Wannier-Stark ladder) as 
well as differences and sums of both quantities. For eFd = huo Q both energy scales coincide 
and the distribution function does not show strong features seen at other high fields: this 
is due to the enhanced resonant tunneling from well to well which (i) leads to enhanced 
current (see Fig. 1) and (ii) prevents accumulation/depletion at certain energies, as seen at 
non-resonant applied fields. This shows that in the high-field region the true distribution 
function cannot be approximated by a heated distribution, neither in k nor in g-space. 

For T = 300 K electron heating is less important (Fig. |]a). Thus the linear regime 
extends to higher fields. Furthermore, the phonon-resonance is hardly visible. Fig. f|(b) 
gives the result for a weakly coupled superlattice with A = 4 meV and an increased impurity 
scattering rate l/r imp = 15 ps _1 . As h/r > k B T, A/2 the calculated velocity-field relations 
deviate significantly both in the low-field mobility and the peak position. Again, the phonon- 
resonance is barely visible here, nor in the calculation for T = 77K where the deviation is 
even larger (not shown here). 

Conclusion: We have performed a fully self-consistent quantum transport calculation for 
SLs which covers the whole range of electric fields. A phonon resonance can be identified 
when the potential drop per period equals the optical phonon energy. The results for the 
drift- velocity and the quantum distribution function have been compared with data obtained 
from the semiclassical Boltzmann equation. In wide-band SLs the Boltzmann equation 
gives reliable results concerning linear response at low fields, electron heating at moderate 
fields, and the onset of negative differential conductivity. In contrast, for high electric fields 
or weakly coupled SLs significant differences appear. In this case the quantum nature of 
transport is important and a semiclassical calculation may be seriously in error. We believe 
that an analysis of the kind presented above can be very useful in checking the quality of 
various approximation schemes. 
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FIG. 1. Drift- velocity versus field for a wide-band SL with A = 20.3 meV, l/ri mp = 3 ps _1 for 
T = 77 K. Full line: Calculation by nonequilibrium Green functions. Dashed line: MC-simulation 
of Boltzmann's transport equation. 
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FIG. 2. Electron distribution versus quasimomentum in the Brillouin zone of the miniband 
for different values of k. (Parameters as Fig. |l|) Full line: NGF calculation. Dashed line: 
MC-simulation of BTE. The dotted line shows the thermal distribution oc exp[— E(q, k)//^^] 
with T e = T in (a) and T e = 140 K in (b) for comparison. 
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FIG. 3. Electron distribution versus quasimomentum k (Parameters as Fig. [l]). (a,b,c): Com- 
parison between NGF calculation (full line) and BTE (dashed line); (d): Results from NGF for 
different fields. The BTE result (not shown) resembles the result from (c) for all three fields. 
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FIG. 4. Drift- velocity versus field for: (a) strongly coupled SL as in Fig. |l] but T = 300K. (b) 
weakly coupled SL with A = 4 meV, l/r imp = 15 ps -1 , T = 300K. Full line: NGF calculation. 
Dashed line: MC-simulation of BTE. 



